Pocillopora Juvenile Fusion
Ariana Huffmyer (R Markdown author)
2019
Setup
Set up workspace, set options, and load required packages.
rm(list=ls(all=TRUE)) knitr::opts_chunk$set(root.dir = "~/R",warning=FALSE, message=FALSE)library(gridExtra)
library(multcomp)
library(emmeans)
library(ggplot2)
library(dplyr)
library(effects)
library(plyr)
library(lattice)
library(glmmTMB)
library(effects)
library(lsmeans)
library(MuMIn)
library(MASS)
library(car)
library(FSA)
library(lmerTest)
library(lme4)
library(blmeco)
library(ggsci)
library(coxme)
library(survival)
library(cowplot)
library(RColorBrewer)
library(tidyverse)
library(partitionBEFsp)
library(factoextra)
library(ggfortify)
library(cowplot)
library(survminer)1. GROW OUT SURVIVORSHIP
Load data
Survivorship is measured as binary on each individual spat in the experiment.
#load data
rearing<-read.csv("Data/GrowOut_Responses.csv", header=TRUE, sep=",", na.strings="NA")
rearing$Quadrant<-as.factor(rearing$Quadrant)
rearing$Diversity.Type<-relevel(rearing$Diversity.Type, ref="Individual")
#rearing$Total.Juveniles<-as.factor(rearing$Total.Juveniles)
#rearing$Fusion.Partners<-as.factor(rearing$Fusion.Partners)
#rearing$Richness<-as.factor(rearing$Richness)Subset the final timepoint for additional analyses.
final<-rearing[which(rearing$Timepoint == "Final"),]
multiple<-rearing[which(rearing$Community == "Multiple"),]
multiple_final<-multiple[which(multiple$Timepoint == "Final"),]Analysis
Analyze with a logistic binomial regression model.
Analyze at the final timepoint.
#with final dataset
survmod2<-glmer(Survivorship ~ Richness + Community + Fusion.Partners + Richness:Fusion.Partners + (1|Parent.Site/Genotype) + (1|Tank) + (1|Slide), data=final, family=binomial)
#final dataset
summary(effect(c("Richness", "Fusion.Partners"), xlevels=4, survmod2))##
## Richness*Fusion.Partners effect
## Fusion.Partners
## Richness 0 1 2 3
## 1 0.4749728 0.9356204 0.9957349 0.9997334
## 2 0.5151037 0.9175013 0.9914845 0.9991803
## 3 0.5550408 0.8948559 0.9830701 0.9974823
## 4 0.5942795 0.8668961 0.9666214 0.9922938
##
## Lower 95 Percent Confidence Limits
## Fusion.Partners
## Richness 0 1 2 3
## 1 0.2857530 0.7943532 0.9407520 0.9836550
## 2 0.3608520 0.8075419 0.9468818 0.9860349
## 3 0.3649140 0.7750230 0.9272927 0.9767284
## 4 0.3218885 0.6699355 0.8147196 0.8865651
##
## Upper 95 Percent Confidence Limits
## Fusion.Partners
## Richness 0 1 2 3
## 1 0.6716629 0.9820395 0.9997088 0.9999957
## 2 0.6665306 0.9671889 0.9986868 0.9999525
## 3 0.7303142 0.9545994 0.9962318 0.9997327
## 4 0.8188353 0.9543351 0.9947841 0.9995289summary(effect("Fusion.Partners", xlevels=4, survmod2))##
## Fusion.Partners effect
## Fusion.Partners
## 0 1 2 3
## 0.5167374 0.9166729 0.9912419 0.9991419
##
## Lower 95 Percent Confidence Limits
## Fusion.Partners
## 0 1 2 3
## 0.3625326 0.8073334 0.9467831 0.9859895
##
## Upper 95 Percent Confidence Limits
## Fusion.Partners
## 0 1 2 3
## 0.6678192 0.9665337 0.9986131 0.9999481summary(effect("Community", xlevels=4, survmod2))##
## Community effect
## Community
## Individual Multiple
## 0.9580409 0.7382699
##
## Lower 95 Percent Confidence Limits
## Community
## Individual Multiple
## 0.8545984 0.5848937
##
## Upper 95 Percent Confidence Limits
## Community
## Individual Multiple
## 0.9888517 0.8495526dispersion_glmer(survmod2) #no overdispersion## [1] 0.9350818plot(residuals(survmod2)) + abline(h=0, lty=2) #Dispersed randomly
## integer(0)Create plot for effect of community type.
eff.surv1 <- predictorEffect("Community", survmod2)
surv.effects1<-plot(eff.surv1,
lines=list(multiline=TRUE,
col=c("black"),
lty=1),
confint=list(style="bars", width=4, col="black"),
lwd=4,
axes=list(y=list(lim=c(0, 1))),
type="response",
ylab=expression(bold("Prob(Survivorship)")),
legend.position="right",
xlab=expression(bold("Number of Juveniles")),
main=""); surv.effects1
#lattice=list(key.args=list(space="right",
#border=FALSE,
#title=expression(bold("Genotypic \nRichness")),
#cex=1,
#cex.title=1)));surv.effects1Create plots for fusion partner effects.
survmod2a<-glmer(Survivorship ~ Richness + Fusion.Partners + Richness:Fusion.Partners + (1|Parent.Site/Genotype) + (1|Tank) + (1|Slide), data=final, family=binomial, subset=c(Community=="Multiple"))
eff.surv2 <- predictorEffect("Fusion.Partners", survmod2a, xlevels=4, rug=TRUE)
surv.effects2<-plot(eff.surv2,
lines=list(multiline=TRUE,
col=c("lightblue","mediumblue", "blue4", "darkblue"),
lty=1),
confint=list(style="bands", alpha=0.1),
lwd=4,
axes=list(y=list(lim=c(0, 1))),
type="response",
ylab=expression(bold("Prob(Survivorship)")),
legend.position="right",
xlab=expression(bold("Fusion Partners")),
main="",
lattice=list(key.args=list(space="right",
border=FALSE,
title=expression(bold("Genotypic \nRichness")),
cex=1,
cex.title=1)));surv.effects2
Generate mean values of survivorship for 1) number of juveniles and 2) within “multiple” juvenile groups, mean values for number of fusion partners.
meansurv <- ddply(final, c("Community"), summarise,
N = length(Survivorship[!is.na(Survivorship)]),
mean = mean(Survivorship, na.rm=TRUE),
sd = sd(Survivorship, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Survivorship, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meansurv## Community N mean sd se max lower upper
## 1 Individual 58 0.7758621 0.4206554 0.05523477 1 0.3552066 1.196518
## 2 Multiple 163 0.6687117 0.4721270 0.03697984 1 0.1965847 1.140839meansurv2 <- ddply(multiple_final, c("Fusion.Partners"), summarise,
N = length(Survivorship[!is.na(Survivorship)]),
mean = mean(Survivorship, na.rm=TRUE),
sd = sd(Survivorship, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Survivorship, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meansurv2## Fusion.Partners N mean sd se max lower upper
## 1 0 75 0.3466667 0.4791133 0.05532324 1 -0.1324466 0.825780
## 2 1 48 0.9583333 0.2019409 0.02914766 1 0.7563924 1.160274
## 3 2 24 0.9166667 0.2823299 0.05763034 1 0.6343368 1.198997
## 4 3 16 0.9375000 0.2500000 0.06250000 1 0.6875000 1.1875002. GROW OUT FUSION
Analysis
Analyze with a logistic binomial regression model.
Final timepoint dataset.
Analyze fusion as a product of genotypic richness within “multiple” juvenile groups.
#with final dataset
fusemod1<-glmer(Fusion ~ Richness + (1|Parent.Site/Genotype) + (1|Tank) + (1|Slide), data=final, family=binomial, subset=c(Community=="Multiple"))
summary(fusemod1)
Anova(fusemod1, type="II") # no effects on rate of fusion by the end point
dispersion_glmer(fusemod1) #no overdispersion
plot(residuals(fusemod1)) + abline(h=0, lty=2)
There is an effect of genotypic diversity on fusion rates. Plot the relationship.
eff.fuse1 <- predictorEffect("Richness", fusemod1, xlevels=4, rug=TRUE)
fuse.effects1<-plot(eff.fuse1,
lines=list(multiline=TRUE,
col=c("black"),
lty=1),
confint=list(style="bands", alpha=0.1),
lwd=4,
#axes=list(y=list(lim=c(0, 1))),
type="response",
ylab=expression(bold("Prob(Successful Fusion)")),
xlab=expression(bold("Richness")),
main="",
lattice=list(key.args=list(space="right",
border=FALSE,
title=expression(bold("Genotypic \nRichness")),
cex=1,
cex.title=1)));fuse.effects1
Generate a summary of proportion fused at the end of the grow-out period.
meanfusion <- ddply(multiple_final, c("Richness"), summarise,
N = length(Fusion[!is.na(Fusion)]),
mean = mean(Fusion, na.rm=TRUE),
sd = sd(Fusion, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Fusion, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meanfusion## Richness N mean sd se max lower upper
## 1 1 57 0.5263158 0.5037454 0.06672270 1 0.02257042 1.0300612
## 2 2 26 0.6153846 0.4961389 0.09730085 1 0.11924568 1.1115236
## 3 3 36 0.3888889 0.4944132 0.08240221 1 -0.10552434 0.8833021
## 4 4 44 0.6363636 0.4866071 0.07335878 1 0.14975654 1.12297073. GROW OUT GROWTH
Analysis
Growth is measured as polyp expansion in each spat by number of spat in each sibling type and number of fused individuals. Growth is measured in final dataset as all spat started with 1 polyp each.
Analyze polyp expansion with the final dataset and present means.
polypmod1<-glmer(Polyps ~ Richness + Community + Fusion.Partners + Richness:Fusion.Partners + (1|Parent.Site/Genotype) + (1|Tank) + (1|Slide), data=final, family = poisson)
summary(polypmod1)## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: poisson ( log )
## Formula:
## Polyps ~ Richness + Community + Fusion.Partners + Richness:Fusion.Partners +
## (1 | Parent.Site/Genotype) + (1 | Tank) + (1 | Slide)
## Data: final
##
## AIC BIC logLik deviance df.resid
## 588.8 616.3 -285.4 570.8 148
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.11373 -0.54677 -0.02775 0.43720 1.40877
##
## Random effects:
## Groups Name Variance Std.Dev.
## Slide (Intercept) 0.000e+00 0.000e+00
## Genotype:Parent.Site (Intercept) 7.465e-02 2.732e-01
## Parent.Site (Intercept) 1.629e-10 1.276e-05
## Tank (Intercept) 9.270e-03 9.628e-02
## Number of obs: 157, groups:
## Slide, 92; Genotype:Parent.Site, 20; Parent.Site, 6; Tank, 5
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.22325 0.14229 8.597 <2e-16 ***
## Richness 0.01432 0.07830 0.183 0.855
## CommunityMultiple -0.10036 0.17493 -0.574 0.566
## Fusion.Partners -0.03007 0.16617 -0.181 0.856
## Richness:Fusion.Partners 0.01653 0.04814 0.343 0.731
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Rchnss CmmntM Fsn.Pr
## Richness -0.545
## CmmntyMltpl 0.073 -0.687
## Fusn.Prtnrs -0.337 0.651 -0.723
## Rchnss:Fs.P 0.408 -0.783 0.669 -0.934
## convergence code: 0
## boundary (singular) fit: see ?isSingularAnova(polypmod1, type="II")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Polyps
## Chisq Df Pr(>Chisq)
## Richness 0.5259 1 0.4683
## Community 0.3291 1 0.5662
## Fusion.Partners 0.1522 1 0.6964
## Richness:Fusion.Partners 0.1179 1 0.7314qqPlot(residuals(polypmod1))
## 813 817
## 144 145meangrowth <- ddply(final, c("Community"), summarise,
N = length(Polyps[!is.na(Polyps)]),
mean = mean(Polyps, na.rm=TRUE),
sd = sd(Polyps, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Polyps, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meangrowth## Community N mean sd se max lower upper
## 1 Individual 45 3.533333 1.778661 0.2651472 8 1.754672 5.311995
## 2 Multiple 112 3.357143 1.558800 0.1472928 7 1.798343 4.915943Plot the model outputs.
eff.polyps1 <- predictorEffect("Fusion.Partners", polypmod1, xlevels=4, rug=TRUE)
polyp.effects1<-plot(eff.polyps1,
lines=list(multiline=TRUE,
col=c("lightblue","mediumblue", "blue4", "darkblue"),
lty=1),
confint=list(style="bands", alpha=0.1),
lwd=4,
#axes=list(y=list(lim=c(0, 1))),
type="response",
ylab=expression(bold("Total Polyp Expansion")),
legend.position="right",
xlab=expression(bold("Fusion Partners")),
main="",
lattice=list(key.args=list(space="right",
border=FALSE,
title=expression(bold("Genotypic \nRichness")),
cex=1,
cex.title=1)));polyp.effects1
Combine figures for Grow-Out Period.
grow_out_figs<-plot_grid(surv.effects1, surv.effects2, labels = c("A", "B"), ncol=2, nrow=1, rel_heights= c(1,1), rel_widths = c(0.8,1), label_size = 20, label_y=1, align="h");grow_out_figs
ggsave(filename="Figures/grow_out_figs.png", plot=grow_out_figs, dpi=500, width=10, height=6, units="in")4. STRESS SURVIVORSHIP
Load data
Survivorship is measured as binary on each individual spat in the experiment.
#load data
stress<-read.csv("Data/Stress_Responses.csv", header=TRUE, sep=",", na.strings="NA")
stress$Diversity.Type<-relevel(stress$Diversity.Type, ref="Individual")
stress$Polyps<-as.numeric(stress$Polyps)Subset the final timepoint, “multiple” juvenile datasets for additional analyses.
final_stress<-stress[which(stress$Timepoint == "End"),]
multiple_stress<-stress[which(stress$Community == "Multiple"),]
multiple_final_stress<-multiple_stress[which(multiple_stress$Timepoint == "End"),]
d26_stress<-stress[which(stress$Timepoint == "S2D11"),]Analysis
Analyze with a logistic binomial regression model.
Check for effect of community type (“individual” or “multiple” juveniles) on survivorship and display means.
#with full time dataset
survmod3b<-glmer(Survivorship ~ Day * Treatment * Community + (1|Parent.Site/Genotype) + (1|Treatment:Tank) + (1|Slide), data=stress, family=binomial, link=logit)
summary(survmod3b)## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula:
## Survivorship ~ Day * Treatment * Community + (1 | Parent.Site/Genotype) +
## (1 | Treatment:Tank) + (1 | Slide)
## Data: stress
##
## AIC BIC logLik deviance df.resid
## 446.3 506.2 -211.2 422.3 1076
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -129.711 0.001 0.026 0.116 6.235
##
## Random effects:
## Groups Name Variance Std.Dev.
## Slide (Intercept) 2.523e+00 1.588500
## Genotype:Parent.Site (Intercept) 5.680e-01 0.753637
## Treatment:Tank (Intercept) 3.662e-06 0.001914
## Parent.Site (Intercept) 8.331e-04 0.028863
## Number of obs: 1088, groups:
## Slide, 90; Genotype:Parent.Site, 20; Treatment:Tank, 9; Parent.Site, 6
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 6.57183 1.66342 3.951 7.79e-05 ***
## Day -0.08040 0.05789 -1.389 0.164885
## TreatmentHigh 0.97861 1.88333 0.520 0.603332
## CommunityMultiple 3.65815 2.44562 1.496 0.134708
## Day:TreatmentHigh -0.28864 0.08140 -3.546 0.000391 ***
## Day:CommunityMultiple -0.18917 0.08695 -2.176 0.029577 *
## TreatmentHigh:CommunityMultiple 3.15953 3.23542 0.977 0.328794
## Day:TreatmentHigh:CommunityMultiple -0.02401 0.12456 -0.193 0.847163
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Day TrtmnH CmmntM Dy:TrH Dy:CmM TrH:CM
## Day -0.832
## TreatmntHgh -0.687 0.656
## CmmntyMltpl -0.571 0.525 0.470
## Dy:TrtmntHg 0.384 -0.624 -0.838 -0.339
## Dy:CmmntyMl 0.494 -0.641 -0.439 -0.918 0.463
## TrtmntHg:CM 0.390 -0.376 -0.579 -0.742 0.487 0.687
## Dy:TrtmH:CM -0.278 0.418 0.545 0.624 -0.627 -0.694 -0.928
## convergence code: 0
## Model failed to converge with max|grad| = 0.0634371 (tol = 0.001, component 1)Anova(survmod3b, type="II")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Survivorship
## Chisq Df Pr(>Chisq)
## Day 60.6800 1 6.715e-15 ***
## Treatment 38.5770 1 5.264e-10 ***
## Community 2.7557 1 0.096909 .
## Day:Treatment 22.1358 1 2.540e-06 ***
## Day:Community 10.3021 1 0.001329 **
## Treatment:Community 4.5938 1 0.032088 *
## Day:Treatment:Community 0.0371 1 0.847163
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1dispersion_glmer(survmod3b) #no overdispersion## [1] 0.558754plot(residuals(survmod3b)) + abline(h=0, lty=2) #Dispersed randomly
## integer(0)meanstresssurv1 <- ddply(final_stress, c("Community", "Treatment"), summarise,
N = length(Survivorship[!is.na(Survivorship)]),
mean = mean(Survivorship, na.rm=TRUE),
sd = sd(Survivorship, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Survivorship, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meanstresssurv1## Community Treatment N mean sd se max lower
## 1 Individual Ambient 21 0.95238095 0.2182179 0.04761905 1 0.7341631
## 2 Individual High 24 0.00000000 0.0000000 0.00000000 0 0.0000000
## 3 Multiple Ambient 46 0.76086957 0.4312660 0.06358670 1 0.3296036
## 4 Multiple High 65 0.04615385 0.2114510 0.02622727 1 -0.1652972
## upper
## 1 1.1705988
## 2 0.0000000
## 3 1.1921355
## 4 0.2576049Analyze with full dataset.
#with full time dataset
survmod3<-glmer(Survivorship ~ Day * Treatment * Fusion.Partners * Richness + (1|Parent.Site/Genotype) + (1|Treatment:Tank) + (1|Slide/Sample), data=stress, family=binomial, subset=c(Community=="Multiple"))
fixed_form<-glm(Survivorship ~ Day * Treatment * Fusion.Partners, data=stress, family=binomial, subset=c(Community=="Multiple"))
library("brglm2"); glm(fixed_form, data = stress, family = binomial, method="detect_separation")## Separation: FALSE
## Existence of maximum likelihood estimates
## (Intercept) Day
## 0 0
## TreatmentHigh Fusion.Partners
## 0 0
## Day:TreatmentHigh Day:Fusion.Partners
## 0 0
## TreatmentHigh:Fusion.Partners Day:TreatmentHigh:Fusion.Partners
## 0 0
## 0: finite value, Inf: infinity, -Inf: -infinitysummary(survmod3)## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula: Survivorship ~ Day * Treatment * Fusion.Partners * Richness +
## (1 | Parent.Site/Genotype) + (1 | Treatment:Tank) + (1 | Slide/Sample)
## Data: stress
## Subset: c(Community == "Multiple")
##
## AIC BIC logLik deviance df.resid
## 266.3 364.0 -112.2 224.3 752
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -9.4233 0.0000 0.0006 0.0429 1.3089
##
## Random effects:
## Groups Name Variance Std.Dev.
## Sample:Slide (Intercept) 6.397e+00 2.5292012
## Slide (Intercept) 1.177e-04 0.0108482
## Genotype:Parent.Site (Intercept) 8.211e-03 0.0906157
## Treatment:Tank (Intercept) 6.207e-04 0.0249143
## Parent.Site (Intercept) 1.662e-07 0.0004077
## Number of obs: 773, groups:
## Sample:Slide, 111; Slide, 45; Genotype:Parent.Site, 20; Treatment:Tank, 8; Parent.Site, 6
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 9.36537 9.77915 0.958 0.3382
## Day -0.36928 0.33310 -1.109 0.2676
## TreatmentHigh -2.26802 10.86374 -0.209 0.8346
## Fusion.Partners 6.29842 7.40902 0.850 0.3953
## Richness 0.91692 7.10758 0.129 0.8974
## Day:TreatmentHigh 0.02190 0.38514 0.057 0.9546
## Day:Fusion.Partners -0.13414 0.24757 -0.542 0.5879
## TreatmentHigh:Fusion.Partners 32.15393 17.20786 1.869 0.0617
## Day:Richness 0.05145 0.23893 0.215 0.8295
## TreatmentHigh:Richness 2.42121 7.54725 0.321 0.7484
## Fusion.Partners:Richness -1.25492 2.99110 -0.420 0.6748
## Day:TreatmentHigh:Fusion.Partners -1.42437 0.69415 -2.052 0.0402
## Day:TreatmentHigh:Richness -0.16342 0.25928 -0.630 0.5285
## Day:Fusion.Partners:Richness 0.01143 0.10035 0.114 0.9093
## TreatmentHigh:Fusion.Partners:Richness -5.63393 5.44568 -1.035 0.3009
## Day:TreatmentHigh:Fusion.Partners:Richness 0.26613 0.21337 1.247 0.2123
##
## (Intercept)
## Day
## TreatmentHigh
## Fusion.Partners
## Richness
## Day:TreatmentHigh
## Day:Fusion.Partners
## TreatmentHigh:Fusion.Partners .
## Day:Richness
## TreatmentHigh:Richness
## Fusion.Partners:Richness
## Day:TreatmentHigh:Fusion.Partners *
## Day:TreatmentHigh:Richness
## Day:Fusion.Partners:Richness
## TreatmentHigh:Fusion.Partners:Richness
## Day:TreatmentHigh:Fusion.Partners:Richness
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## convergence code: 0
## unable to evaluate scaled gradient
## Model failed to converge: degenerate Hessian with 2 negative eigenvalues
## failure to converge in 10000 evaluationsAnova(survmod3, type="II")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Survivorship
## Chisq Df Pr(>Chisq)
## Day 52.1650 1 5.103e-13 ***
## Treatment 36.7599 1 1.336e-09 ***
## Fusion.Partners 0.0832 1 0.77299
## Richness 1.4063 1 0.23567
## Day:Treatment 5.9263 1 0.01492 *
## Day:Fusion.Partners 5.5057 1 0.01895 *
## Treatment:Fusion.Partners 1.0661 1 0.30183
## Day:Richness 0.1740 1 0.67658
## Treatment:Richness 0.0571 1 0.81121
## Fusion.Partners:Richness 1.1707 1 0.27927
## Day:Treatment:Fusion.Partners 5.1440 1 0.02333 *
## Day:Treatment:Richness 0.0168 1 0.89683
## Day:Fusion.Partners:Richness 0.6301 1 0.42733
## Treatment:Fusion.Partners:Richness 1.3224 1 0.25016
## Day:Treatment:Fusion.Partners:Richness 1.5557 1 0.21230
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1dispersion_glmer(survmod3) #no overdispersion## [1] 0.4350411plot(residuals(survmod3)) + abline(h=0, lty=2) #Dispersed randomly
## integer(0)meanstresssurv2 <- ddply(final_stress, c("Fusion.Partners", "Treatment"), summarise,
N = length(Survivorship[!is.na(Survivorship)]),
mean = mean(Survivorship, na.rm=TRUE),
sd = sd(Survivorship, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Survivorship, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meanstresssurv2## Fusion.Partners Treatment N mean sd se max lower
## 1 0 Ambient 30 0.83333333 0.3790490 0.06920457 1 0.4542843
## 2 0 High 50 0.04000000 0.1979487 0.02799417 1 -0.1579487
## 3 1 Ambient 16 0.81250000 0.4031129 0.10077822 1 0.4093871
## 4 1 High 15 0.06666667 0.2581989 0.06666667 1 -0.1915322
## 5 2 Ambient 9 0.77777778 0.4409586 0.14698618 1 0.3368192
## 6 2 High 16 0.00000000 0.0000000 0.00000000 0 0.0000000
## 7 3 Ambient 12 0.83333333 0.3892495 0.11236664 1 0.4440839
## 8 3 High 8 0.00000000 0.0000000 0.00000000 0 0.0000000
## upper
## 1 1.2123824
## 2 0.2379487
## 3 1.2156129
## 4 0.3248656
## 5 1.2187363
## 6 0.0000000
## 7 1.2225828
## 8 0.0000000Create plot for effect of number of juveniles (“individual” or “multiple”).
eff.surv3b <- predictorEffect("Community", survmod3b)
surv.effects3b<-plot(eff.surv3b,
lines=list(multiline=TRUE,
col=c("blue", "red"),
lty=1),
confint=list(style="bars", width=4, col="black"),
lwd=4,
#axes=list(y=list(lim=c(0, 1))),
type="rescale",
ylab=expression(bold("Prob(Survivorship)")),
legend.position="none",
xlab=expression(bold("Number of Juveniles")),
main="",
lattice=list(key.args=list(space="none",
border=FALSE,
title=expression(bold("Temperature")),
cex=1,
cex.title=1)));surv.effects3b
Plot effect of number of juveniles over time.
stress$group<-paste(stress$Treatment, "-", stress$Community)
survmod3b2<-glmer(Survivorship ~ group * Day + (1|Parent.Site/Genotype) + (1|Treatment:Tank) + (1|Slide), data=stress, family=binomial, link=logit)
eff.surv3b2 <- predictorEffect(c("Day"), survmod3b2)
surv.effects3b2<-plot(eff.surv3b2,
lines=list(multiline=TRUE,
col=c("blue", "blue", "red", "red"),
lty=c(1, 2, 1, 2)),
confint=list(style="bands", alpha=0.1),
lwd=4,
axes=list(y=list(lim=c(0, 1.1))),
type="response",
ylab=expression(bold("Prob(Survivorship)")),
legend.position="top",
xlab=expression(bold("Day")),
main="",
lattice=list(key.args=list(space="top",
border=FALSE,
title=expression(bold("Temperature - Number of Juveniles")),
cex=1,
cex.title=1)));surv.effects3b2
Examine curve estimates to identify 50% survivorship.
summary(effect(c("group", "Day"), xlevels=50, survmod3b2))##
## group*Day effect
## Day
## group 0 0.653 1.31 1.96 2.61
## Ambient - Individual 0.9986417 0.9985680 0.9984898 0.9984082 0.9983223
## Ambient - Multiple 0.9999603 0.9999528 0.9999437 0.9999331 0.9999204
## High - Individual 0.9994462 0.9992965 0.9991051 0.9988646 0.9985596
## High - Multiple 0.9999994 0.9999992 0.9999988 0.9999982 0.9999974
## Day
## group 3.27 3.92 4.57 5.22 5.88
## Ambient - Individual 0.9982303 0.9981347 0.9980340 0.9979279 0.9978142
## Ambient - Multiple 0.9999051 0.9998871 0.9998658 0.9998403 0.9998096
## High - Individual 0.9981660 0.9976737 0.9970497 0.9962588 0.9952396
## High - Multiple 0.9999961 0.9999943 0.9999917 0.9999879 0.9999823
## Day
## group 6.53 7.18 7.84 8.49 9.14
## Ambient - Individual 0.9976963 0.9975720 0.9974389 0.9973007 0.9971551
## Ambient - Multiple 0.9997736 0.9997308 0.9996790 0.9996182 0.9995460
## High - Individual 0.9939665 0.9923555 0.9902834 0.9877011 0.9844434
## High - Multiple 0.9999741 0.9999622 0.9999445 0.9999190 0.9998817
## Day
## group 9.8 10.4 11.1 11.8 12.4
## Ambient - Individual 0.9969992 0.9968502 0.9966669 0.9964730 0.9962979
## Ambient - Multiple 0.9994586 0.9993647 0.9992345 0.9990775 0.9989175
## High - Individual 0.9802693 0.9755340 0.9686015 0.9597856 0.9503831
## High - Multiple 0.9998262 0.9997536 0.9996297 0.9994434 0.9992109
## Day
## group 13.1 13.7 14.4 15 15.7
## Ambient - Individual 0.9960827 0.9958883 0.9956493 0.9954335 0.9951682
## Ambient - Multiple 0.9986957 0.9984696 0.9981561 0.9978368 0.9973940
## High - Individual 0.9367858 0.9224404 0.9019753 0.8807360 0.8510437
## High - Multiple 0.9988143 0.9983194 0.9974758 0.9964241 0.9946345
## Day
## group 16.3 17 17.6 18.3 18.9
## Ambient - Individual 0.9949286 0.9946341 0.9943683 0.9940415 0.9937464
## Ambient - Multiple 0.9969431 0.9963180 0.9956818 0.9948000 0.9939032
## High - Individual 0.8209593 0.7801006 0.7400643 0.6877657 0.6387043
## High - Multiple 0.9924079 0.9886311 0.9839535 0.9760725 0.9664033
## Day
## group 19.6 20.2 20.9 21.6 22.2
## Ambient - Individual 0.9933837 0.9930563 0.9926540 0.9922285 0.9918443
## Ambient - Multiple 0.9926609 0.9913983 0.9896510 0.9875533 0.9854248
## High - Individual 0.5776507 0.5232799 0.4592354 0.3965098 0.3452513
## High - Multiple 0.9503378 0.9310043 0.8997671 0.8565663 0.8080989
## Day
## group 22.9 23.5 24.2 24.8 25.5
## Ambient - Individual 0.9913724 0.9909463 0.9904229 0.9899504 0.9893700
## Ambient - Multiple 0.9824856 0.9795085 0.9754059 0.9712606 0.9655648
## High - Individual 0.2897524 0.2466542 0.2021122 0.1689489 0.1359078
## High - Multiple 0.7369400 0.6639104 0.5678757 0.4809678 0.3813680
## Day
## group 26.1 26.8 27.4 28.1 28.7
## Ambient - Individual 0.9888462 0.98820274 0.98762213 0.98690902 0.98626562
## Ambient - Multiple 0.9598291 0.95197972 0.94411219 0.93340482 0.92274101
## High - Individual 0.1120816 0.08897129 0.07268148 0.05717206 0.04640781
## High - Multiple 0.3029899 0.22431753 0.16937894 0.11945350 0.08730690
## Day
## group 29.4 30 30.7 31.3 32
## Ambient - Individual 0.98547551 0.98476275 0.98388759 0.98309824 0.98212922
## Ambient - Multiple 0.90833673 0.89411415 0.87509511 0.85652939 0.83202769
## High - Individual 0.03628552 0.02933141 0.02284452 0.01841715 0.01430846
## High - Multiple 0.05983017 0.04294646 0.02898720 0.02061633 0.01381049
##
## Lower 95 Percent Confidence Limits
## Day
## group 0 0.653 1.31 1.96 2.61
## Ambient - Individual 0.9646655 0.9649766 0.9652683 0.9655355 0.9657805
## Ambient - Multiple 0.9979830 0.9977864 0.9975689 0.9973324 0.9970723
## High - Individual 0.9913895 0.9898742 0.9880798 0.9859913 0.9835369
## High - Multiple 0.9999769 0.9999692 0.9999588 0.9999452 0.9999270
## Day
## group 3.27 3.92 4.57 5.22 5.88
## Ambient - Individual 0.9660059 0.9662040 0.9663773 0.9665247 0.9666467
## Ambient - Multiple 0.9967818 0.9964668 0.9961203 0.9957389 0.9953123
## High - Individual 0.9806048 0.9772087 0.9732211 0.9685412 0.9629614
## High - Multiple 0.9999023 0.9998698 0.9998265 0.9997688 0.9996904
## Day
## group 6.53 7.18 7.84 8.49 9.14
## Ambient - Individual 0.9667381 0.9667996 0.9668295 0.9668250 0.9667847
## Ambient - Multiple 0.9948491 0.9943388 0.9937671 0.9931458 0.9924601
## High - Individual 0.9565144 0.9489683 0.9400006 0.9296820 0.9176650
## High - Multiple 0.9995872 0.9994494 0.9992622 0.9990153 0.9986853
## Day
## group 9.8 10.4 11.1 11.8 12.4
## Ambient - Individual 0.9667047 0.9665953 0.9664205 0.9661908 0.9659467
## Ambient - Multiple 0.9916908 0.9909205 0.9899270 0.9888192 0.9877678
## High - Individual 0.9034731 0.8885987 0.8685655 0.8453005 0.8225148
## High - Multiple 0.9982364 0.9976957 0.9968508 0.9956941 0.9943680
## Day
## group 13.1 13.7 14.4 15 15.7
## Ambient - Individual 0.9656022 0.9652515 0.9647720 0.9642954 0.9636556
## Ambient - Multiple 0.9864083 0.9851157 0.9834407 0.9818446 0.9797712
## High - Individual 0.7923188 0.7631391 0.7250822 0.6889804 0.6429070
## High - Multiple 0.9922933 0.9899136 0.9861903 0.9819210 0.9752488
## Day
## group 16.3 17 17.6 18.3 18.9
## Ambient - Individual 0.9630289 0.9621975 0.9613909 0.9603290 0.9593052
## Ambient - Multiple 0.9777900 0.9752087 0.9727341 0.9694981 0.9663838
## High - Individual 0.6002734 0.5474025 0.5000354 0.4433987 0.3946586
## High - Multiple 0.9676146 0.9557286 0.9422083 0.9213483 0.8979280
## Day
## group 19.6 20.2 20.9 21.6 22.2
## Ambient - Individual 0.9579646 0.9566780 0.9549996 0.9531103 0.9513053
## Ambient - Multiple 0.9622932 0.9583378 0.9531159 0.9471326 0.9412965
## High - Individual 0.3388896 0.2931041 0.2432415 0.1982310 0.1640523
## High - Multiple 0.8624812 0.8237211 0.7671948 0.6975931 0.6279100
## Day
## group 22.9 23.5 24.2 24.8 25.5
## Ambient - Individual 0.9489604 0.9467251 0.94382720 0.94107013 0.93750221
## Ambient - Multiple 0.9335201 0.9258872 0.91565244 0.90554622 0.89192097
## High - Individual 0.1295765 0.1046011 0.08047226 0.06366912 0.04798639
## High - Multiple 0.5372850 0.4551979 0.36018275 0.28439130 0.20743679
## Day
## group 26.1 26.8 27.4 28.1 28.7
## Ambient - Individual 0.93411387 0.92973700 0.92558810 0.92023915 0.91517905
## Ambient - Multiple 0.87840706 0.86013421 0.84199477 0.81751490 0.79334071
## High - Individual 0.03739302 0.02775739 0.02139171 0.01570692 0.01200962
## High - Multiple 0.15343341 0.10461494 0.07369029 0.04796865 0.03274520
## Day
## group 29.4 30 30.7 31.3
## Ambient - Individual 0.908669500 0.902525769 0.894642221 0.887222149
## Ambient - Multiple 0.761024773 0.729555825 0.688298487 0.649103361
## High - Individual 0.008749816 0.006652628 0.004819988 0.003649879
## High - Multiple 0.020718341 0.013881109 0.008636886 0.005723046
## Day
## group 32
## Ambient - Individual 0.877729690
## Ambient - Multiple 0.599268815
## High - Individual 0.002633742
## High - Multiple 0.003525362
##
## Upper 95 Percent Confidence Limits
## Day
## group 0 0.653 1.31 1.96 2.61
## Ambient - Individual 0.9999495 0.9999433 0.9999364 0.9999288 0.9999203
## Ambient - Multiple 0.9999992 0.9999990 0.9999987 0.9999983 0.9999978
## High - Individual 0.9999647 0.9999516 0.9999335 0.9999091 0.9998757
## High - Multiple 1.0000000 1.0000000 1.0000000 0.9999999 0.9999999
## Day
## group 3.27 3.92 4.57 5.22 5.88
## Ambient - Individual 0.9999107 0.9999002 0.9998885 0.9998755 0.9998609
## Ambient - Multiple 0.9999972 0.9999964 0.9999954 0.9999940 0.9999923
## High - Individual 0.9998294 0.9997669 0.9996819 0.9995660 0.9994055
## High - Multiple 0.9999998 0.9999998 0.9999996 0.9999994 0.9999990
## Day
## group 6.53 7.18 7.84 8.49 9.14
## Ambient - Individual 0.9998451 0.9998275 0.9998079 0.9997866 0.9997631
## Ambient - Multiple 0.9999901 0.9999873 0.9999836 0.9999789 0.9999728
## High - Individual 0.9991902 0.9988977 0.9984940 0.9979542 0.9972245
## High - Multiple 0.9999984 0.9999974 0.9999958 0.9999933 0.9999894
## Day
## group 9.8 10.4 11.1 11.8 12.4
## Ambient - Individual 0.9997371 0.9997112 0.9996782 0.9996421 0.9996085
## Ambient - Multiple 0.9999650 0.9999559 0.9999423 0.9999246 0.9999052
## High - Individual 0.9962224 0.9950079 0.9931037 0.9904985 0.9875264
## High - Multiple 0.9999829 0.9999737 0.9999566 0.9999283 0.9998899
## Day
## group 13.1 13.7 14.4 15 15.7
## Ambient - Individual 0.9995660 0.9995267 0.9994773 0.9994319 0.9993753
## Ambient - Multiple 0.9998762 0.9998445 0.9997974 0.9997459 0.9996695
## High - Individual 0.9829246 0.9777298 0.9697905 0.9609650 0.9477285
## High - Multiple 0.9998186 0.9997219 0.9995429 0.9993010 0.9988547
## Day
## group 16.3 17 17.6 18.3 18.9
## Ambient - Individual 0.9993237 0.9992598 0.9992019 0.9991310 0.9990673
## Ambient - Multiple 0.9995863 0.9994630 0.9993294 0.9991323 0.9989194
## High - Individual 0.9333370 0.9123221 0.8901685 0.8589709 0.8273928
## High - Multiple 0.9982544 0.9971533 0.9956827 0.9930096 0.9894800
## Day
## group 19.6 20.2 20.9 21.6 22.2
## Ambient - Individual 0.9989901 0.9989215 0.9988391 0.9987546 0.9986808
## Ambient - Multiple 0.9986070 0.9982714 0.9977819 0.9971622 0.9965044
## High - Individual 0.7849112 0.7439744 0.6917146 0.6358343 0.5862351
## High - Multiple 0.9831613 0.9749784 0.9607114 0.9392466 0.9131048
## Day
## group 22.9 23.5 24.2 24.8 25.5
## Ambient - Individual 0.9985938 0.9985188 0.9984314 0.9983570 0.9982714
## Ambient - Multiple 0.9955574 0.9945621 0.9931457 0.9916758 0.9896129
## High - Individual 0.5278548 0.4785237 0.4230312 0.3780295 0.3292136
## High - Multiple 0.8711125 0.8236434 0.7541648 0.6836186 0.5921726
## Day
## group 26.1 26.8 27.4 28.1 28.7
## Ambient - Individual 0.9981994 0.9981177 0.9980500 0.9979741 0.9979120
## Ambient - Multiple 0.9875043 0.9845934 0.9816688 0.9777043 0.9737934
## High - Individual 0.2908738 0.2504116 0.2193783 0.1872750 0.1630685
## High - Multiple 0.5104289 0.4171719 0.3432751 0.2675319 0.2127824
## Day
## group 29.4 30 30.7 31.3 32
## Ambient - Individual 0.9978434 0.9977881 0.99772791 0.99768010 0.99762885
## Ambient - Multiple 0.9685890 0.9635462 0.95694978 0.95065978 0.94255216
## High - Individual 0.1383791 0.1199838 0.10140441 0.08767469 0.07389981
## High - Multiple 0.1606633 0.1251475 0.09279867 0.07148059 0.05252057Create plot for effect of temperature on survival in multiple juvenile groups across number of fusion partners.
stress$group<-paste(stress$Treatment, "-", stress$Fusion.Partners)
survmod3a<-glmer(Survivorship ~ Day * group + (1|Parent.Site/Genotype) + (1|Treatment:Tank) + (1|Slide), data=stress, family=binomial, link=logit)
summary(effect(c("Day", "group"), xlevels=50, survmod3a))##
## Day*group effect
## group
## Day Ambient - 0 Ambient - 1 Ambient - 2 Ambient - 3 High - 0 High - 1
## 0 0.9996572 0.9999989 0.9999970 0.9999981 0.99982061 1.000000e+00
## 0.653 0.9996213 0.9999986 0.9999963 0.9999976 0.99977069 1.000000e+00
## 1.31 0.9995814 0.9999983 0.9999954 0.9999971 0.99970645 1.000000e+00
## 1.96 0.9995378 0.9999978 0.9999942 0.9999964 0.99962520 1.000000e+00
## 2.61 0.9994896 0.9999972 0.9999928 0.9999955 0.99952147 1.000000e+00
## 3.27 0.9994355 0.9999964 0.9999911 0.9999944 0.99938676 1.000000e+00
## 3.92 0.9993767 0.9999955 0.9999889 0.9999931 0.99921710 1.000000e+00
## 4.57 0.9993118 0.9999943 0.9999862 0.9999914 0.99900055 1.000000e+00
## 5.22 0.9992400 0.9999927 0.9999828 0.9999894 0.99872418 1.000000e+00
## 5.88 0.9991596 0.9999907 0.9999785 0.9999868 0.99836538 1.000000e+00
## 6.53 0.9990720 0.9999883 0.9999733 0.9999837 0.99791373 1.000000e+00
## 7.18 0.9989753 0.9999851 0.9999668 0.9999798 0.99733763 1.000000e+00
## 7.84 0.9988669 0.9999811 0.9999585 0.9999749 0.99659022 1.000000e+00
## 8.49 0.9987489 0.9999760 0.9999484 0.9999690 0.99565023 1.000000e+00
## 9.14 0.9986186 0.9999696 0.9999359 0.9999616 0.99445256 1.000000e+00
## 9.8 0.9984724 0.9999614 0.9999200 0.9999522 0.99290099 9.999999e-01
## 10.4 0.9983262 0.9999519 0.9999021 0.9999418 0.99112022 9.999999e-01
## 11.1 0.9981378 0.9999379 0.9998762 0.9999268 0.98847721 9.999997e-01
## 11.8 0.9979283 0.9999199 0.9998435 0.9999078 0.98505939 9.999993e-01
## 12.4 0.9977301 0.9999003 0.9998086 0.9998877 0.98134852 9.999985e-01
## 13.1 0.9974748 0.9998713 0.9997579 0.9998587 0.97586784 9.999964e-01
## 13.7 0.9972334 0.9998399 0.9997040 0.9998278 0.96994367 9.999923e-01
## 14.4 0.9969224 0.9997934 0.9996256 0.9997833 0.96124408 9.999818e-01
## 15 0.9966283 0.9997428 0.9995422 0.9997360 0.95190686 9.999616e-01
## 15.7 0.9962496 0.9996681 0.9994210 0.9996677 0.93831863 9.999087e-01
## 16.3 0.9958914 0.9995870 0.9992920 0.9995953 0.92389542 9.998079e-01
## 17 0.9954303 0.9994670 0.9991047 0.9994906 0.90319768 9.995427e-01
## 17.6 0.9949943 0.9993367 0.9989053 0.9993796 0.88159829 9.990384e-01
## 18.3 0.9944331 0.9991441 0.9986159 0.9992192 0.85124931 9.977130e-01
## 18.9 0.9939025 0.9989350 0.9983077 0.9990490 0.82036409 9.952009e-01
## 19.6 0.9932198 0.9986258 0.9978607 0.9988032 0.77826709 9.886469e-01
## 20.2 0.9925744 0.9982904 0.9973848 0.9985426 0.73691184 9.764115e-01
## 20.9 0.9917442 0.9977945 0.9966948 0.9981661 0.68281940 9.456000e-01
## 21.6 0.9908219 0.9971551 0.9958234 0.9976925 0.62329066 8.795086e-01
## 22.2 0.9899507 0.9964619 0.9948967 0.9971907 0.56903668 7.762713e-01
## 22.9 0.9888303 0.9954380 0.9935543 0.9964661 0.50367541 5.930040e-01
## 23.5 0.9877724 0.9943289 0.9921282 0.9956987 0.44746580 4.091888e-01
## 24.2 0.9864125 0.9926922 0.9900649 0.9945915 0.38363796 2.253093e-01
## 24.8 0.9851290 0.9909216 0.9878761 0.9934198 0.33186713 1.214567e-01
## 25.5 0.9834801 0.9883132 0.9847154 0.9917308 0.27628347 5.486888e-02
## 26.1 0.9819246 0.9854971 0.9813702 0.9899455 0.23351129 2.685464e-02
## 26.8 0.9799276 0.9813593 0.9765538 0.9873761 0.18972323 1.145549e-02
## 27.4 0.9780451 0.9769067 0.9714741 0.9846651 0.15743652 5.478218e-03
## 28.1 0.9756302 0.9703918 0.9641925 0.9807725 0.12557683 2.307795e-03
## 28.7 0.9733556 0.9634171 0.9565539 0.9766769 0.10282118 1.098325e-03
## 29.4 0.9704405 0.9532786 0.9456767 0.9708167 0.08095186 4.615125e-04
## 30 0.9676976 0.9425114 0.9343572 0.9646768 0.06567531 2.194302e-04
## 30.7 0.9641861 0.9270180 0.9183968 0.9559376 0.05125529 9.215677e-05
## 31.3 0.9608862 0.9107653 0.9019821 0.9468387 0.04133080 4.380825e-05
## 32 0.9566675 0.8877354 0.8791680 0.9339881 0.03207250 1.839681e-05
## group
## Day High - 2 High - 3
## 0 1.000000e+00 1.000000e+00
## 0.653 1.000000e+00 1.000000e+00
## 1.31 1.000000e+00 1.000000e+00
## 1.96 1.000000e+00 1.000000e+00
## 2.61 1.000000e+00 1.000000e+00
## 3.27 1.000000e+00 1.000000e+00
## 3.92 1.000000e+00 1.000000e+00
## 4.57 1.000000e+00 1.000000e+00
## 5.22 1.000000e+00 1.000000e+00
## 5.88 1.000000e+00 1.000000e+00
## 6.53 1.000000e+00 1.000000e+00
## 7.18 1.000000e+00 1.000000e+00
## 7.84 1.000000e+00 1.000000e+00
## 8.49 1.000000e+00 1.000000e+00
## 9.14 1.000000e+00 1.000000e+00
## 9.8 1.000000e+00 1.000000e+00
## 10.4 1.000000e+00 1.000000e+00
## 11.1 1.000000e+00 1.000000e+00
## 11.8 1.000000e+00 1.000000e+00
## 12.4 1.000000e+00 1.000000e+00
## 13.1 1.000000e+00 1.000000e+00
## 13.7 1.000000e+00 1.000000e+00
## 14.4 1.000000e+00 9.999999e-01
## 15 1.000000e+00 9.999999e-01
## 15.7 1.000000e+00 9.999996e-01
## 16.3 1.000000e+00 9.999989e-01
## 17 9.999999e-01 9.999967e-01
## 17.6 9.999996e-01 9.999915e-01
## 18.3 9.999982e-01 9.999744e-01
## 18.9 9.999937e-01 9.999345e-01
## 19.6 9.999719e-01 9.998034e-01
## 20.2 9.998995e-01 9.994961e-01
## 20.9 9.995553e-01 9.984897e-01
## 21.6 9.980336e-01 9.954827e-01
## 22.2 9.929953e-01 9.884978e-01
## 22.9 9.697086e-01 9.662694e-01
## 23.5 8.994109e-01 9.178410e-01
## 24.2 6.687835e-01 7.883074e-01
## 24.8 3.606032e-01 5.922039e-01
## 25.5 1.129702e-01 3.261762e-01
## 26.1 3.435032e-02 1.587978e-01
## 26.8 7.968983e-03 5.919969e-02
## 27.4 2.238669e-03 2.395144e-02
## 28.1 5.064186e-04 8.113347e-03
## 28.7 1.414991e-04 3.179750e-03
## 29.4 3.195713e-05 1.062166e-03
## 30 8.926134e-06 4.144881e-04
## 30.7 2.015732e-06 1.382006e-04
## 31.3 5.630138e-07 5.389950e-05
## 32 1.271410e-07 1.796711e-05
##
## Lower 95 Percent Confidence Limits
## group
## Day Ambient - 0 Ambient - 1 Ambient - 2 Ambient - 3 High - 0
## 0 0.9912518 0.9545168 0.7325815 0.8462365 0.997701705
## 0.653 0.9908699 0.9539437 0.7364089 0.8484868 0.997242915
## 1.31 0.9904657 0.9533518 0.7401424 0.8506719 0.996688404
## 1.96 0.9900450 0.9527502 0.7437172 0.8527544 0.996029605
## 2.61 0.9896022 0.9521316 0.7471703 0.8547566 0.995239134
## 3.27 0.9891284 0.9514850 0.7505482 0.8567057 0.994274598
## 3.92 0.9886362 0.9508288 0.7537444 0.8585408 0.993133072
## 4.57 0.9881170 0.9501518 0.7568066 0.8602898 0.991763107
## 5.22 0.9875687 0.9494527 0.7597298 0.8619502 0.990118988
## 5.88 0.9869801 0.9487184 0.7625502 0.8635426 0.988112714
## 6.53 0.9863666 0.9479692 0.7651760 0.8650155 0.985738645
## 7.18 0.9857169 0.9471922 0.7676444 0.8663899 0.982890630
## 7.84 0.9850174 0.9463723 0.7699816 0.8676802 0.979417641
## 8.49 0.9842859 0.9455317 0.7721083 0.8688425 0.975312202
## 9.14 0.9835088 0.9446553 0.7740511 0.8698915 0.970394334
## 9.8 0.9826691 0.9437253 0.7758245 0.8708339 0.964409325
## 10.4 0.9818578 0.9428414 0.7772512 0.8715769 0.957941763
## 11.1 0.9808483 0.9417593 0.7786770 0.8722975 0.948930606
## 11.8 0.9797648 0.9406163 0.7798260 0.8728488 0.938048523
## 12.4 0.9787716 0.9395827 0.7805725 0.8731755 0.926966414
## 13.1 0.9775303 0.9383067 0.7811410 0.8733718 0.911640875
## 13.7 0.9763890 0.9371458 0.7813462 0.8733674 0.896137579
## 14.4 0.9749581 0.9357031 0.7812235 0.8731403 0.874877034
## 15 0.9736381 0.9343812 0.7807770 0.8727359 0.853587459
## 15.7 0.9719774 0.9327256 0.7798132 0.8719913 0.824756426
## 16.3 0.9704400 0.9311960 0.7785645 0.8710918 0.796316308
## 17 0.9684988 0.9292627 0.7765530 0.8696980 0.758497888
## 17.6 0.9666950 0.9274589 0.7742936 0.8681694 0.721983661
## 18.3 0.9644087 0.9251547 0.7709464 0.8659395 0.674651080
## 18.9 0.9622761 0.9229797 0.7673835 0.8635894 0.630273632
## 19.6 0.9595626 0.9201658 0.7622936 0.8602527 0.574670885
## 20.2 0.9570216 0.9174732 0.7570106 0.8568006 0.524491444
## 20.9 0.9537761 0.9139373 0.7495969 0.8519589 0.464245411
## 21.6 0.9501908 0.9098839 0.7406070 0.8460764 0.403748426
## 22.2 0.9468142 0.9059037 0.7314288 0.8400456 0.353041394
## 22.9 0.9424772 0.9005273 0.7186987 0.8316244 0.296738636
## 23.5 0.9383789 0.8951399 0.7057623 0.8229870 0.251991951
## 24.2 0.9330984 0.8877021 0.6878924 0.8109048 0.204863053
## 24.8 0.9280945 0.8800777 0.6698116 0.7984833 0.169291314
## 25.5 0.9216313 0.8692986 0.6449788 0.7810695 0.133604370
## 26.1 0.9154940 0.8579844 0.6200568 0.7631454 0.107867452
## 26.8 0.9075543 0.8416163 0.5862431 0.7380417 0.083081857
## 27.4 0.9000072 0.8240793 0.5528906 0.7123116 0.065853898
## 28.1 0.8902394 0.7982953 0.5087392 0.6766100 0.049782817
## 28.7 0.8809566 0.7704196 0.4665948 0.6405886 0.038919461
## 29.4 0.8689542 0.7295108 0.4131631 0.5918255 0.029020218
## 30 0.8575680 0.6859797 0.3648143 0.5442980 0.022462141
## 30.7 0.8428857 0.6242950 0.3073540 0.4829192 0.016584592
## 31.3 0.8290071 0.5622163 0.2590662 0.4265387 0.012745536
## 32 0.8111942 0.4810510 0.2061730 0.3588126 0.009344534
## group
## Day High - 1 High - 2 High - 3
## 0 9.999827e-01 1.000000e+00 9.998425e-01
## 0.653 9.999758e-01 1.000000e+00 9.998034e-01
## 1.31 9.999662e-01 1.000000e+00 9.997543e-01
## 1.96 9.999530e-01 1.000000e+00 9.996936e-01
## 2.61 9.999346e-01 1.000000e+00 9.996177e-01
## 3.27 9.999084e-01 1.000000e+00 9.995213e-01
## 3.92 9.998725e-01 1.000000e+00 9.994025e-01
## 4.57 9.998223e-01 1.000000e+00 9.992540e-01
## 5.22 9.997525e-01 1.000000e+00 9.990682e-01
## 5.88 9.996533e-01 9.999999e-01 9.988318e-01
## 6.53 9.995167e-01 9.999999e-01 9.985398e-01
## 7.18 9.993260e-01 9.999998e-01 9.981740e-01
## 7.84 9.990552e-01 9.999996e-01 9.977077e-01
## 8.49 9.986818e-01 9.999993e-01 9.971307e-01
## 9.14 9.981604e-01 9.999987e-01 9.964066e-01
## 9.8 9.974187e-01 9.999977e-01 9.954813e-01
## 10.4 9.964869e-01 9.999960e-01 9.944316e-01
## 11.1 9.949653e-01 9.999924e-01 9.928895e-01
## 11.8 9.927832e-01 9.999856e-01 9.909120e-01
## 12.4 9.901732e-01 9.999750e-01 9.887758e-01
## 13.1 9.859144e-01 9.999524e-01 9.856261e-01
## 13.7 9.808288e-01 9.999172e-01 9.822144e-01
## 14.4 9.725550e-01 9.998419e-01 9.771691e-01
## 15 9.627213e-01 9.997243e-01 9.716861e-01
## 15.7 9.468436e-01 9.994717e-01 9.635476e-01
## 16.3 9.281818e-01 9.990757e-01 9.546667e-01
## 17 8.985543e-01 9.982206e-01 9.414220e-01
## 17.6 8.645462e-01 9.968729e-01 9.268921e-01
## 18.3 8.123640e-01 9.939432e-01 9.050880e-01
## 18.9 7.551412e-01 9.892942e-01 8.809994e-01
## 19.6 6.726272e-01 9.791248e-01 8.445556e-01
## 20.2 5.889855e-01 9.629366e-01 8.039375e-01
## 20.9 4.797506e-01 9.277187e-01 7.419533e-01
## 21.6 3.648355e-01 8.606862e-01 6.591930e-01
## 22.2 2.690919e-01 7.616605e-01 5.673259e-01
## 22.9 1.698139e-01 5.817398e-01 4.331900e-01
## 23.5 1.020533e-01 3.838546e-01 3.015108e-01
## 24.2 4.799674e-02 1.695646e-01 1.563212e-01
## 24.8 2.183644e-02 5.899768e-02 6.913025e-02
## 25.5 7.583508e-03 1.161293e-02 1.995638e-02
## 26.1 2.808945e-03 2.266247e-03 5.660537e-03
## 26.8 8.250626e-04 2.859700e-04 1.117317e-03
## 27.4 2.780836e-04 4.474087e-05 2.552580e-04
## 28.1 7.598612e-05 4.866239e-06 4.274697e-05
## 28.7 2.457363e-05 7.058022e-07 8.901133e-06
## 29.4 6.494592e-06 7.260488e-08 1.384993e-06
## 30 2.058156e-06 1.020469e-08 2.760031e-07
## 30.7 5.346023e-07 1.023491e-09 4.139013e-08
## 31.3 1.675413e-07 1.416432e-10 8.058400e-09
## 32 4.308669e-08 1.401818e-11 1.183879e-09
##
## Upper 95 Percent Confidence Limits
## group
## Day Ambient - 0 Ambient - 1 Ambient - 2 Ambient - 3 High - 0 High - 1
## 0 0.9999867 1.0000000 1.0000000 1.0000000 0.9999860 1.000000000
## 0.653 0.9999844 1.0000000 1.0000000 1.0000000 0.9999810 1.000000000
## 1.31 0.9999818 1.0000000 1.0000000 1.0000000 0.9999740 1.000000000
## 1.96 0.9999787 1.0000000 1.0000000 1.0000000 0.9999647 1.000000000
## 2.61 0.9999752 1.0000000 1.0000000 1.0000000 0.9999521 1.000000000
## 3.27 0.9999710 1.0000000 1.0000000 1.0000000 0.9999346 1.000000000
## 3.92 0.9999662 1.0000000 1.0000000 1.0000000 0.9999112 1.000000000
## 4.57 0.9999606 1.0000000 1.0000000 1.0000000 0.9998795 1.000000000
## 5.22 0.9999541 1.0000000 1.0000000 1.0000000 0.9998365 1.000000000
## 5.88 0.9999464 1.0000000 1.0000000 1.0000000 0.9997772 1.000000000
## 6.53 0.9999376 1.0000000 1.0000000 1.0000000 0.9996980 1.000000000
## 7.18 0.9999274 1.0000000 1.0000000 1.0000000 0.9995908 1.000000000
## 7.84 0.9999154 1.0000000 1.0000000 1.0000000 0.9994433 1.000000000
## 8.49 0.9999017 1.0000000 1.0000000 1.0000000 0.9992466 1.000000000
## 9.14 0.9998859 1.0000000 1.0000000 1.0000000 0.9989811 1.000000000
## 9.8 0.9998673 1.0000000 1.0000000 1.0000000 0.9986167 1.000000000
## 10.4 0.9998479 1.0000000 1.0000000 1.0000000 0.9981751 1.000000000
## 11.1 0.9998218 0.9999999 0.9999999 1.0000000 0.9974814 1.000000000
## 11.8 0.9997914 0.9999999 0.9999999 0.9999999 0.9965288 1.000000000
## 12.4 0.9997614 0.9999998 0.9999999 0.9999999 0.9954361 1.000000000
## 13.1 0.9997213 0.9999997 0.9999998 0.9999999 0.9937302 0.999999999
## 13.7 0.9996818 0.9999996 0.9999997 0.9999998 0.9917830 0.999999997
## 14.4 0.9996291 0.9999994 0.9999995 0.9999997 0.9887614 0.999999988
## 15 0.9995775 0.9999991 0.9999993 0.9999995 0.9853367 0.999999962
## 15.7 0.9995087 0.9999985 0.9999988 0.9999992 0.9800682 0.999999851
## 16.3 0.9994416 0.9999977 0.9999982 0.9999989 0.9741575 0.999999523
## 17 0.9993525 0.9999963 0.9999972 0.9999983 0.9651785 0.999998146
## 17.6 0.9992659 0.9999944 0.9999959 0.9999975 0.9552544 0.999994086
## 18.3 0.9991516 0.9999909 0.9999935 0.9999961 0.9404515 0.999977252
## 18.9 0.9990409 0.9999864 0.9999905 0.9999943 0.9244386 0.999928290
## 19.6 0.9988954 0.9999782 0.9999853 0.9999912 0.9011669 0.999729128
## 20.2 0.9987553 0.9999674 0.9999786 0.9999873 0.8767399 0.999164361
## 20.9 0.9985722 0.9999481 0.9999671 0.9999806 0.8424778 0.996957269
## 21.6 0.9983658 0.9999178 0.9999498 0.9999706 0.8016993 0.989334387
## 22.2 0.9981689 0.9998786 0.9999283 0.9999583 0.7616126 0.970326274
## 22.9 0.9979138 0.9998099 0.9998925 0.9999379 0.7093606 0.912115385
## 23.5 0.9976719 0.9997224 0.9998490 0.9999132 0.6606486 0.808450940
## 24.2 0.9973606 0.9995718 0.9997781 0.9998732 0.6005848 0.626552626
## 24.8 0.9970674 0.9993844 0.9996945 0.9998262 0.5476434 0.461248455
## 25.5 0.9966928 0.9990708 0.9995625 0.9997520 0.4858822 0.306063530
## 26.1 0.9963424 0.9986933 0.9994122 0.9996677 0.4342658 0.212812347
## 26.8 0.9958979 0.9980865 0.9991839 0.9995397 0.3769704 0.139877979
## 27.4 0.9954850 0.9973892 0.9989349 0.9993998 0.3312229 0.098353582
## 28.1 0.9949650 0.9963290 0.9985738 0.9991965 0.2824639 0.065778486
## 28.7 0.9944854 0.9951847 0.9981987 0.9989846 0.2449067 0.046889950
## 29.4 0.9938854 0.9935632 0.9976822 0.9986915 0.2060904 0.031782362
## 30 0.9933358 0.9919383 0.9971732 0.9984011 0.1769724 0.022869494
## 30.7 0.9926526 0.9898059 0.9965089 0.9980197 0.1475328 0.015640733
## 31.3 0.9920307 0.9878220 0.9958880 0.9976607 0.1258535 0.011326114
## 32 0.9912622 0.9853919 0.9951180 0.9972124 0.1042620 0.007793989
## group
## Day High - 2 High - 3
## 0 1.000000000 1.0000000
## 0.653 1.000000000 1.0000000
## 1.31 1.000000000 1.0000000
## 1.96 1.000000000 1.0000000
## 2.61 1.000000000 1.0000000
## 3.27 1.000000000 1.0000000
## 3.92 1.000000000 1.0000000
## 4.57 1.000000000 1.0000000
## 5.22 1.000000000 1.0000000
## 5.88 1.000000000 1.0000000
## 6.53 1.000000000 1.0000000
## 7.18 1.000000000 1.0000000
## 7.84 1.000000000 1.0000000
## 8.49 1.000000000 1.0000000
## 9.14 1.000000000 1.0000000
## 9.8 1.000000000 1.0000000
## 10.4 1.000000000 1.0000000
## 11.1 1.000000000 1.0000000
## 11.8 1.000000000 1.0000000
## 12.4 1.000000000 1.0000000
## 13.1 1.000000000 1.0000000
## 13.7 1.000000000 1.0000000
## 14.4 1.000000000 1.0000000
## 15 1.000000000 1.0000000
## 15.7 1.000000000 1.0000000
## 16.3 1.000000000 1.0000000
## 17 1.000000000 1.0000000
## 17.6 1.000000000 1.0000000
## 18.3 0.999999999 1.0000000
## 18.9 0.999999996 1.0000000
## 19.6 0.999999963 0.9999998
## 20.2 0.999999738 0.9999990
## 20.9 0.999997459 0.9999934
## 21.6 0.999976017 0.9999602
## 22.2 0.999841004 0.9998225
## 22.9 0.998644654 0.9990696
## 23.5 0.992267885 0.9965532
## 24.2 0.952306511 0.9868145
## 24.8 0.835339624 0.9659833
## 25.5 0.579920206 0.9200471
## 26.1 0.357778717 0.8622565
## 26.8 0.184063272 0.7797273
## 27.4 0.101134071 0.7022452
## 28.1 0.050111396 0.6101571
## 28.7 0.027592771 0.5333971
## 29.4 0.013871742 0.4494359
## 30 0.007747416 0.3838474
## 30.7 0.003954239 0.3158071
## 31.3 0.002232915 0.2650040
## 32 0.001151806 0.2142610eff.surv3 <- predictorEffect(c("Day"), survmod3a)
surv.effects3<-plot(eff.surv3,
lines=list(multiline=TRUE,
col=c("blue", "blue", "blue", "blue", "red", "red", "red", "red"),
lty=c(1, 2, 4, 3, 1, 2, 4, 3)),
confint=list(style="bands", alpha=0.1),
lwd=4,
axes=list(y=list(lim=c(0, 1.1))),
type="response",
ylab=expression(bold("Prob(Survivorship)")),
legend.position="top",
xlab=expression(bold("Day")),
main="",
lattice=list(key.args=list(space="top",
border=FALSE,
title=expression(bold("Temperature - Fusion Partners")),
cex=1,
cex.title=1)));surv.effects3
Display mean survivorship within multiple juvenile groups at the end of the experiment.
meanstresssurv2 <- ddply(multiple_final_stress, c("Fusion.Partners", "Treatment"), summarise,
N = length(Survivorship[!is.na(Survivorship)]),
mean = mean(Survivorship, na.rm=TRUE),
sd = sd(Survivorship, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Survivorship, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meanstresssurv2## Fusion.Partners Treatment N mean sd se max lower
## 1 0 Ambient 9 0.55555556 0.5270463 0.17568209 1 0.02850928
## 2 0 High 26 0.07692308 0.2717465 0.05329387 1 -0.19482341
## 3 1 Ambient 16 0.81250000 0.4031129 0.10077822 1 0.40938711
## 4 1 High 15 0.06666667 0.2581989 0.06666667 1 -0.19153222
## 5 2 Ambient 9 0.77777778 0.4409586 0.14698618 1 0.33681923
## 6 2 High 16 0.00000000 0.0000000 0.00000000 0 0.00000000
## 7 3 Ambient 12 0.83333333 0.3892495 0.11236664 1 0.44408386
## 8 3 High 8 0.00000000 0.0000000 0.00000000 0 0.00000000
## upper
## 1 1.0826018
## 2 0.3486696
## 3 1.2156129
## 4 0.3248656
## 5 1.2187363
## 6 0.0000000
## 7 1.2225828
## 8 0.00000005. STRESS FUSION
Analysis
Analyze with a logistic binomial regression model.
Analyze fusion success over the course of the stress test.
#with full
fusemod2<-glmer(Fusion ~ Day * Treatment * Richness + (1|Parent.Site/Genotype) + (1|Treatment:Tank) + (1|Slide/Sample), data=stress, family=binomial, link=logit, subset=c(Community=="Multiple"))
summary(fusemod2)
Anova(fusemod2, type="II")
plot(allEffects(fusemod2))
dispersion_glmer(fusemod2) #no overdispersion
plot(residuals(fusemod2)) + abline(h=0, lty=2)
There is an effect of treatment and day and richness. Plot this relationship.
stress$group<-paste(stress$Treatment, "-", stress$Richness)
fusemod2a<-glmer(Fusion ~ Day * group + (1|Parent.Site/Genotype) + (1|Treatment:Tank) + (1|Slide/Sample), data=stress, family=binomial, link=logit, subset=c(Community=="Multiple"))
eff.fuse2 <- predictorEffect("Day", xlevels=4, fusemod2a, rug=TRUE)
fuse.effects2<-plot(eff.fuse2,
lines=list(multiline=TRUE,
col=c("blue", "blue", "blue", "blue", "red", "red", "red", "red"),
lty=c(1, 2, 4, 3, 1, 2, 4, 3)),
confint=list(style="none"),
lwd=4,
axes=list(y=list(lim=c(0, 1))),
type="response",
ylab=expression(bold("Prob(Successful Fusion)")),
legend.position="right",
xlab=expression(bold("Time (Days)")),
main="",
lattice=list(key.args=list(space="none",
border=FALSE,
title=expression(bold("Temperature - \nGenotypic Richness")),
cex=1,
cex.title=1)));fuse.effects2
fusemod2b<-glmer(Fusion ~ Day * Treatment + (1|Parent.Site/Genotype) + (1|Treatment:Tank) + (1|Slide/Sample), data=stress, family=binomial, link=logit, subset=c(Community=="Multiple"))
Anova(fusemod2b, type="II") ## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Fusion
## Chisq Df Pr(>Chisq)
## Day 32.6383 1 1.11e-08 ***
## Treatment 0.2578 1 0.6117
## Day:Treatment 0.0725 1 0.7877
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1eff.fuse2b <- predictorEffect("Day", xlevels=4, fusemod2b, rug=TRUE)
fuse.effects2b<-plot(eff.fuse2b,
lines=list(multiline=TRUE,
col=c("blue", "red"),
lty=c(1)),
confint=list(style="bands"),
lwd=4,
axes=list(y=list(lim=c(0.8, 1))),
type="response",
ylab=expression(bold("Prob(Successful Fusion)")),
legend.position="right",
xlab=expression(bold("Time (Days)")),
main="",
lattice=list(key.args=list(space="none",
border=FALSE,
title=expression(bold("Temperature")),
cex=1,
cex.title=1)));fuse.effects2b
Analyze fusion as a product of treatment and genotypic richness at the end of the experiment.
#with final dataset
fusemod3<-glmer(Fusion ~ Treatment + (1|Treatment:Tank) + (1|Parent.Site/Genotype) + (1|Slide), data=stress, link=logit, family=binomial, subset=c(Community=="Multiple" & Timepoint == "End"))
summary(fusemod3)## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula:
## Fusion ~ Treatment + (1 | Treatment:Tank) + (1 | Parent.Site/Genotype) +
## (1 | Slide)
## Data: stress
## Subset: c(Community == "Multiple" & Timepoint == "End")
##
## AIC BIC logLik deviance df.resid
## 117.0 133.2 -52.5 105.0 105
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.43276 -0.20590 0.04881 0.15110 1.32248
##
## Random effects:
## Groups Name Variance Std.Dev.
## Slide (Intercept) 40.4085 6.3568
## Genotype:Parent.Site (Intercept) 0.8975 0.9474
## Treatment:Tank (Intercept) 0.0000 0.0000
## Parent.Site (Intercept) 0.0000 0.0000
## Number of obs: 111, groups:
## Slide, 45; Genotype:Parent.Site, 20; Treatment:Tank, 8; Parent.Site, 6
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.532 2.268 2.439 0.0147 *
## TreatmentHigh -4.525 2.568 -1.762 0.0781 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## TreatmntHgh -0.658
## convergence code: 0
## boundary (singular) fit: see ?isSingularAnova(fusemod3, type="II") # no effects on rate of fusion by the end point## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Fusion
## Chisq Df Pr(>Chisq)
## Treatment 3.1045 1 0.07808 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(allEffects(fusemod3))
dispersion_glmer(fusemod3) #no overdispersion## [1] 0.6023656plot(residuals(fusemod3)) + abline(h=0, lty=2)
## integer(0)meanstressfuse <- ddply(multiple_final_stress, c("Treatment", "Richness", "Timepoint"), summarise,
N = length(Fusion[!is.na(Fusion)]),
mean = mean(Fusion, na.rm=TRUE),
sd = sd(Fusion, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Fusion, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meanstressfuse## Treatment Richness Timepoint N mean sd se max lower
## 1 Ambient 1 End 17 0.7058824 0.4696682 0.1139113 1 0.23621413
## 2 Ambient 2 End 14 0.7142857 0.4688072 0.1252940 1 0.24547848
## 3 Ambient 3 End 3 1.0000000 0.0000000 0.0000000 1 1.00000000
## 4 Ambient 4 End 12 1.0000000 0.0000000 0.0000000 1 1.00000000
## 5 High 1 End 25 0.6000000 0.5000000 0.1000000 1 0.10000000
## 6 High 2 End 6 0.6666667 0.5163978 0.2108185 1 0.15026889
## 7 High 3 End 18 0.5000000 0.5144958 0.1212678 1 -0.01449576
## 8 High 4 End 16 0.6875000 0.4787136 0.1196784 1 0.20878645
## upper
## 1 1.175551
## 2 1.183093
## 3 1.000000
## 4 1.000000
## 5 1.100000
## 6 1.183064
## 7 1.014496
## 8 1.166214eff.fuse3 <- predictorEffect("Treatment", fusemod3, xlevels=4, rug=TRUE)
fuse.effects3<-plot(eff.fuse3,
lines=list(multiline=TRUE,
col=c("black"),
lty=c(1,2,4,3)),
confint=list(style="bars"),
lwd=4,
#axes=list(y=list(lim=c(0, 1))),
type="response",
ylab=expression(bold("Prob(Successful Fusion)")),
legend.position="right",
xlab=expression(bold("Temperature Treatment")),
main="",
lattice=list(key.args=list(space="right",
border=FALSE,
title=expression(bold("Genotypic \nRichness")),
cex=1,
cex.title=1)));fuse.effects3
No effects at the end of the experiment, we will report full analysis over time during the stress test.
6. STRESS GROWTH
Analyze with a poisson regression model.
Analyze growth during the course of the stress test.
Test for effect of temperature:
#with full dataset
#polypmod2<-glmer(Polyps ~ Treatment + Community + Community:Treatment + (1|Parent.Site/Genotype) + (1|Treatment:Tank) + (1|Slide), data=stress, family = poisson, subset=c(Timepoint=="S2D11"))
#summary(polypmod2)
#Anova(polypmod2, type="II")
#qqPlot(residuals(polypmod2))
polypmod2a<-glmer(Polyps^2 ~ Day * Treatment * Community + (1|Parent.Site/Genotype) + (1|Treatment:Tank) + (1|Slide/Sample), data=stress, family = poisson)
summary(polypmod2a)
Anova(polypmod2a, type="II")
qqPlot(residuals(polypmod2a))
Plot difference in temperature:
stress$group<-paste(stress$Treatment, "-", stress$Community)
polypmod2b<-glmer(Polyps^2 ~ Day * group + (1|Parent.Site/Genotype) + (1|Treatment:Tank) + (1|Slide/Sample), data=stress, family = poisson)
eff.surv5 <- predictorEffect("Day", xlevels=4, polypmod2b)
polyp.effects2<-plot(eff.surv5,
lines=list(multiline=TRUE,
col=c("blue", "blue", "red", "red"),
lty=c(1,2, 1, 2)),
confint=list(style="bands", width=4),
lwd=4,
#axes=list(y=list(lim=c(0, 1))),
type="response",
ylab=expression(bold("(Total Polyp Expansion)"^2)),
legend.position="right",
xlab=expression(bold("Time (Days)")),
main="",
lattice=list(key.args=list(space="top",
border=FALSE,
title=expression(bold("Treatment - Number of Juveniles")),
cex=1,
cex.title=1)));polyp.effects2
Present means of growth by number of juveniles. Sample size was too low due to mortality at the end of the experiment to calculate growth within multiple juvenile groups, so we pool these together for analysis here.
meanstressgrowth <- ddply(stress, c("Community", "Treatment", "Timepoint"), summarise,
N = length(Polyps[!is.na(Polyps)]),
mean = mean(Polyps, na.rm=TRUE),
sd = sd(Polyps, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Polyps, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meanstressgrowthCombine figures for Stress Period
stress_figs<-plot_grid(fuse.effects2, surv.effects3b2, surv.effects3, polyp.effects2, labels = c("A", "B", "C", "D"), ncol=4, nrow=1, rel_heights= c(1,1,1,1), rel_widths = c(0.6,1,1,1), label_size = 20, label_y=1, align="h");stress_figs
ggsave(filename="Figures/stress_figs.png", plot=stress_figs, dpi=500, width=20, height=6, units="in")7. TEMPERATURE DATA
Load in temperature data from Apex files. Data were recorded every 15 minutes and probes were calibrated to a digital thermometer weekly.
juvtemps<-read.csv("Data/ApexTemps.csv", header=T, na.strings="NA")
#convert Date and Time to date/time format
juvtemps$Date <- as.POSIXct(juvtemps$Date, format="%m/%d/%y %H:%M")
#convert to long format
juvtemps<-juvtemps %>% gather(Tank, Temperature, Tank17:Tank24)
juvtemps$Tank<-as.factor(juvtemps$Tank)
dailymax<-juvtemps[which(juvtemps$TimeDay == "930"),]
dailymax<-juvtemps[which(juvtemps$Phase == "Stress"),]
dailymax$Treatment<-ifelse(dailymax$Tank %in% c("Tank17"), "Ambient",
ifelse(dailymax$Tank %in% c("Tank18"), "Ambient",
ifelse(dailymax$Tank %in% c("Tank19"), "High",
ifelse(dailymax$Tank %in% c("Tank20"), "High",
ifelse(dailymax$Tank %in% c("Tank21"), "Ambient",
ifelse(dailymax$Tank %in% c("Tank22"), "High",
ifelse(dailymax$Tank %in% c("Tank23"), "Ambient",
ifelse(dailymax$Tank %in% c("Tank24"), "High",NA))))))))
maxtemps <- ddply(dailymax, c("Treatment"), summarise,
N = length(Temperature[!is.na(Temperature)]),
mean = mean(Temperature, na.rm=TRUE),
sd = sd(Temperature, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Temperature, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
juvtemps$Period<-ifelse(juvtemps$Phase %in% c("Stress"), "Stress",
ifelse(juvtemps$Phase %in% c("Stress2"), "Stress",
ifelse(juvtemps$Phase %in% c("Acclimation"), "Stress",
ifelse(juvtemps$Phase %in% c("Hurricane"), "Hurricane",
ifelse(juvtemps$Phase %in% c("Ambient"), "Grow-Out",NA)))))Plot timeseries of temperature data with tank temperature in colors.
#subset temperature measurements
juvtempsHIGH<-subset(juvtemps, Period==c("Stress"), c(Date, TimeDay, Phase, Tank, Temperature, Period))
juvtempsHIGH$Treatment<-ifelse(juvtempsHIGH$Tank %in% c("Tank17"), "Ambient",
ifelse(juvtempsHIGH$Tank %in% c("Tank18"), "Ambient",
ifelse(juvtempsHIGH$Tank %in% c("Tank19"), "High",
ifelse(juvtempsHIGH$Tank %in% c("Tank20"), "High",
ifelse(juvtempsHIGH$Tank %in% c("Tank21"), "Ambient",
ifelse(juvtempsHIGH$Tank %in% c("Tank22"), "High",
ifelse(juvtempsHIGH$Tank %in% c("Tank23"), "Ambient",
ifelse(juvtempsHIGH$Tank %in% c("Tank24"), "High",NA))))))))
juvtempsHURR<-subset(juvtemps, Period==c("Hurricane"), c(Date, TimeDay, Phase, Tank, Temperature, Period))
juvtempsHURR$Treatment<-ifelse(juvtempsHURR$Tank %in% c("Tank17"), "Ambient",
ifelse(juvtempsHURR$Tank %in% c("Tank18"), "Ambient",
ifelse(juvtempsHURR$Tank %in% c("Tank19"), "High",
ifelse(juvtempsHURR$Tank %in% c("Tank20"), "High",
ifelse(juvtempsHURR$Tank %in% c("Tank21"), "Ambient",
ifelse(juvtempsHURR$Tank %in% c("Tank22"), "High",
ifelse(juvtempsHURR$Tank %in% c("Tank23"), "Ambient",
ifelse(juvtempsHURR$Tank %in% c("Tank24"), "High",NA))))))))
juvtempsGROW<-subset(juvtemps, Period==c("Grow-Out"), c(Date, TimeDay, Phase, Tank, Temperature, Period))
juvtempsGROW$Treatment<-ifelse(juvtempsGROW$Tank %in% c("Tank17"), "Ambient",
ifelse(juvtempsGROW$Tank %in% c("Tank18"), "Ambient",
ifelse(juvtempsGROW$Tank %in% c("Tank19"), "Ambient",
ifelse(juvtempsGROW$Tank %in% c("Tank20"), "Ambient",
ifelse(juvtempsGROW$Tank %in% c("Tank21"), "Ambient",
ifelse(juvtempsGROW$Tank %in% c("Tank22"), "Ambient",
ifelse(juvtempsGROW$Tank %in% c("Tank23"), "Ambient",
ifelse(juvtempsGROW$Tank %in% c("Tank24"), "Ambient",NA))))))))
juvtemps<-rbind(juvtempsHIGH, juvtempsGROW)
juvtemps<-rbind(juvtemps, juvtempsHURR)
TimeSeries1j<-ggplot(juvtemps, aes(x = Date, y = Temperature)) +
geom_line(aes(color = Treatment), size = 1) +
scale_color_manual(values=c("blue", "red")) +
ylab("Temperature (°C)")+
xlab("Date")+
ylim(26,32)+
geom_vline(xintercept = as.POSIXct(as.Date(c("2018-07-31 10:00:00"))), linetype="solid",
color = "black", size=1)+
geom_vline(xintercept = as.POSIXct(as.Date(c("2018-08-16 00:00:00"))), linetype="dashed",
color = "black", size=1)+
geom_vline(xintercept = as.POSIXct(as.Date(c("2018-09-20 16:00:00"))), linetype="dotted",
color = "black", size=1)+
theme_classic()+
scale_x_datetime(limits = as.POSIXct(as.Date(c("2018-07-30 09:00:00", "2018-09-20 21:00:00")))) +
theme(legend.position="none")+
theme(text = element_text(size = 18, color="black"))+
theme(axis.text = element_text(size=16, color="black"));TimeSeries1j
Plot by time of day by summarizing the mean temp at each time of day for each tank. Display the mean temperatures of each treatment during the Stress test phase (excluding the pause in treatment due to the hurricane system shutdown).
#subset juv temps for the stress period
juvtempsSTRESS<-juvtemps
juvtempsSTRESS$Period<-ifelse(juvtempsSTRESS$Phase %in% c("Stress"), "Stress",
ifelse(juvtempsSTRESS$Phase %in% c("Stress2"), "Stress",
ifelse(juvtempsSTRESS$Phase %in% c("Acclimation"), "Acclimation",
ifelse(juvtempsSTRESS$Phase %in% c("Hurricane"), "Hurricane",
ifelse(juvtempsSTRESS$Phase %in% c("Ambient"), "Grow-Out",NA)))))
juvtempsSTRESS
juvtempsSTRESS<-juvtempsSTRESS[which(juvtempsSTRESS$Phase == "Stress"),]
QuarterlyJ <- ddply(juvtempsSTRESS, c("TimeDay", "Treatment"), summarise,
N = length(Temperature[!is.na(Temperature)]),
mean = mean(Temperature, na.rm=TRUE),
sd = sd(Temperature, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Temperature, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
MeansJtemp <- ddply(juvtemps, c("Period", "Treatment"), summarise,
N = length(Temperature[!is.na(Temperature)]),
mean = mean(Temperature, na.rm=TRUE),
sd = sd(Temperature, na.rm=TRUE),
se = sd / sqrt(N)
)
MeansJtemp
MeansJtemp2 <- ddply(juvtemps, c("Treatment"), summarise,
N = length(Temperature[!is.na(Temperature)]),
mean = mean(Temperature, na.rm=TRUE),
max = max(Temperature, na.rm=TRUE),
sd = sd(Temperature, na.rm=TRUE),
se = sd / sqrt(N)
)
MeansJtemp2Generate final figure.
JuvTempsDay<-ggplot(QuarterlyJ, aes(x = TimeDay, y = mean, color=Treatment, fill=Treatment)) +
geom_line(color="black", size = 1) +
#facet_wrap(~Period) +
scale_x_continuous(breaks=c(0, 720, 1440), labels=c("0:00", "12:00", "24:00")) +
scale_color_manual(values = c("blue", "red")) +
scale_fill_manual(values = c("blue", "red")) +
theme_classic() +
geom_ribbon(aes(ymin=QuarterlyJ$lower, ymax=QuarterlyJ$upper), linetype=2, alpha=0.4, color=NA) +
ylab("Temperature (°C)") +
xlab("Time of Day") +
ylim(26,32)+
theme(text = element_text(size = 18, color="black"))+
theme(panel.margin = unit(2, "lines"))+
theme(axis.text = element_text(size=16, color="black")); JuvTempsDay
Combine figures for temperature.
temp_figs<-plot_grid(TimeSeries1j, JuvTempsDay, labels = c("A", "B"), ncol=2, nrow=1, rel_heights= c(1,1), rel_widths = c(1,0.8), label_size = 20, label_y=1, align="h");temp_figs
ggsave(filename="Figures/temp_figs.png", plot=temp_figs, dpi=500, width=14, height=6, units="in")